FigureS5L

Load Data

dataset <- read.delim("raw_data/FigureS5L.txt", stringsAsFactors = FALSE)

dataset$genotype <- gsub(" ", "", dataset$genotype )
dataset$genotype <- factor(dataset$genotype)
dataset$Experiment <- factor(rep(paste0("exp", 1:(length(dataset$genotype)/length(levels(dataset$genotype)))), each=length(unique(dataset$genotype))))

dataset$UID <- factor(paste(dataset$Experiment, dataset$genotype))

# wide format
kable(dataset, row.names = F)

genotype	NT	bleomycin_5uM	bleomycin_10uM	bleomycin_50uM	Experiment	UID
WT	1160	810	567	151	exp1	exp1 WT
YFP-ALC1	1005	792	585	192	exp1	exp1 YFP-ALC1
WT	1056	838	378	258	exp2	exp2 WT
YFP-ALC1	2132	1742	983	399	exp2	exp2 YFP-ALC1
WT	1604	1036	558	184	exp3	exp3 WT
YFP-ALC1	1714	1325	985	311	exp3	exp3 YFP-ALC1

library(reshape2)
# reshape to long format
dataset <- melt(dataset, variable.name = "Treatment", value.name = "Counts")

dataset$genotype <- relevel(dataset$genotype, ref = "WT")
dataset$UID <- relevel(dataset$UID, ref = "exp1 WT")

dataset$Bleomycin <- gsub("NT","1",dataset$Treatment)
dataset$Bleomycin <- gsub("bleomycin_|uM","",dataset$Bleomycin)
dataset$Bleomycin <- log10(as.integer(dataset$Bleomycin))




dataset$Offset <- NA
for(uid in levels(dataset$UID)){
        dataset$Offset[dataset$UID == uid] <- mean(dataset$Counts[dataset$UID == uid])
}

dataset$NormCounts <- dataset$Counts / dataset$Offset



dataset$Offset2 <- NA
for(gid in levels(dataset$genotype)){
        dataset$Offset2[dataset$genotype == gid] <- mean(dataset$NormCounts[dataset$genotype == gid & dataset$Bleomycin == 0])
}

dataset$NormCounts2 <- dataset$NormCounts / dataset$Offset2



# long format
kable(dataset, row.names = F)

genotype	Experiment	UID	Treatment	Counts	Bleomycin	Offset	NormCounts	Offset2	NormCounts2
WT	exp1	exp1 WT	NT	1160	0.00000	672.00	1.7261905	1.764286	0.9784074
YFP-ALC1	exp1	exp1 YFP-ALC1	NT	1005	0.00000	643.50	1.5617716	1.588615	0.9831029
WT	exp2	exp2 WT	NT	1056	0.00000	632.50	1.6695652	1.764286	0.9463121
YFP-ALC1	exp2	exp2 YFP-ALC1	NT	2132	0.00000	1314.00	1.6225266	1.588615	1.0213469
WT	exp3	exp3 WT	NT	1604	0.00000	845.50	1.8971023	1.764286	1.0752805
YFP-ALC1	exp3	exp3 YFP-ALC1	NT	1714	0.00000	1083.75	1.5815456	1.588615	0.9955502
WT	exp1	exp1 WT	bleomycin_5uM	810	0.69897	672.00	1.2053571	1.764286	0.6831983
YFP-ALC1	exp1	exp1 YFP-ALC1	bleomycin_5uM	792	0.69897	643.50	1.2307692	1.588615	0.7747438
WT	exp2	exp2 WT	bleomycin_5uM	838	0.69897	632.50	1.3249012	1.764286	0.7509560
YFP-ALC1	exp2	exp2 YFP-ALC1	bleomycin_5uM	1742	0.69897	1314.00	1.3257230	1.588615	0.8345152
WT	exp3	exp3 WT	bleomycin_5uM	1036	0.69897	845.50	1.2253105	1.764286	0.6945078
YFP-ALC1	exp3	exp3 YFP-ALC1	bleomycin_5uM	1325	0.69897	1083.75	1.2226067	1.588615	0.7696056
WT	exp1	exp1 WT	bleomycin_10uM	567	1.00000	672.00	0.8437500	1.764286	0.4782388
YFP-ALC1	exp1	exp1 YFP-ALC1	bleomycin_10uM	585	1.00000	643.50	0.9090909	1.588615	0.5722539
WT	exp2	exp2 WT	bleomycin_10uM	378	1.00000	632.50	0.5976285	1.764286	0.3387367
YFP-ALC1	exp2	exp2 YFP-ALC1	bleomycin_10uM	983	1.00000	1314.00	0.7480974	1.588615	0.4709118
WT	exp3	exp3 WT	bleomycin_10uM	558	1.00000	845.50	0.6599645	1.764286	0.3740689
YFP-ALC1	exp3	exp3 YFP-ALC1	bleomycin_10uM	985	1.00000	1083.75	0.9088812	1.588615	0.5721219
WT	exp1	exp1 WT	bleomycin_50uM	151	1.69897	672.00	0.2247024	1.764286	0.1273617
YFP-ALC1	exp1	exp1 YFP-ALC1	bleomycin_50uM	192	1.69897	643.50	0.2983683	1.588615	0.1878167
WT	exp2	exp2 WT	bleomycin_50uM	258	1.69897	632.50	0.4079051	1.764286	0.2312013
YFP-ALC1	exp2	exp2 YFP-ALC1	bleomycin_50uM	399	1.69897	1314.00	0.3036530	1.588615	0.1911433
WT	exp3	exp3 WT	bleomycin_50uM	184	1.69897	845.50	0.2176227	1.764286	0.1233489
YFP-ALC1	exp3	exp3 YFP-ALC1	bleomycin_50uM	311	1.69897	1083.75	0.2869666	1.588615	0.1806395

Plot Data

library(ggplot2)

# raw data
ggplot(dataset, aes(x=Bleomycin, y=Counts)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE, aes(colour=genotype)) +
        geom_point(aes(colour=genotype, shape=Experiment), size=2) +        
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_shape_manual(values=15:20) +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts Linear
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Linear
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ x, se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts Quadratic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)")+
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Quadratic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,2), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts Cubic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)")+
        scale_color_manual(values=c("#000000","#808000"))

# NormCounts2 Cubic
ggplot(dataset, aes(x=Bleomycin, y=NormCounts2, color=genotype)) + 
        theme_bw() +
        theme(panel.grid=element_blank(), text = element_text(size=14)) +
        geom_point(aes(colour=genotype), size=2) +        
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=FALSE) +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        scale_color_manual(values=c("#000000","#808000"))

library(Cairo)

cairo_pdf("FigureS5L.pdf", width = 5, height = 4, family = "Arial")

ggplot(dataset, aes(x=Bleomycin, y=NormCounts2)) + 
        theme_bw() +
        theme(panel.grid.major=element_blank(), panel.grid.minor=element_blank(), 
              axis.line = element_line(colour = "black"), text = element_text(size=14),
              panel.border = element_blank(), panel.background = element_blank()) +
        geom_point(aes(colour = genotype)) +
        geom_smooth(method=lm, formula = y ~ poly(x,3), se=TRUE, aes(colour = genotype), fill='#CCCCCC') +
        #facet_grid(. ~ genotype) +
        xlab(label = "Bleomycin (log10 µM)") +
        ylab(label = "Normalized Counts") +
        scale_color_manual(values=c("#000000","#808000"))

dev.off()

## quartz_off_screen 
##                 2

Models

library(MASS)
library(DHARMa)
library(lme4)
library(lmerTest)
library(bbmle)

Linear formula

fit1 <- lm(Counts ~ Experiment + Bleomycin*genotype, data = dataset)
print(summary(fit1))

## 
## Call:
## lm(formula = Counts ~ Experiment + Bleomycin * genotype, data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -476.82 -107.86  -15.27  130.19  500.52 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  1064.5      141.8   7.509 5.97e-07 ***
## Experimentexp2                315.5      123.7   2.550   0.0201 *  
## Experimentexp3                306.9      123.7   2.481   0.0232 *  
## Bleomycin                    -653.7      117.1  -5.583 2.68e-05 ***
## genotypeYFP-ALC1              417.3      173.2   2.410   0.0269 *  
## Bleomycin:genotypeYFP-ALC1   -141.6      165.6  -0.855   0.4038    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 247.4 on 18 degrees of freedom
## Multiple R-squared:  0.8398, Adjusted R-squared:  0.7953 
## F-statistic: 18.88 on 5 and 18 DF,  p-value: 1.355e-06

cat("AIC: ", AIC(fit1))

## AIC:  339.7383

simres <- simulateResiduals(fittedModel = fit1)
plot(simres)

fit2 <- lm(NormCounts ~ Bleomycin*genotype, data = dataset)
print(summary(fit2))

## 
## Call:
## lm(formula = NormCounts ~ Bleomycin * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.26678 -0.05879 -0.01889  0.09365  0.20857 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 1.76526    0.06366  27.729  < 2e-16 ***
## Bleomycin                  -0.90086    0.06087 -14.799 3.08e-12 ***
## genotypeYFP-ALC1           -0.10406    0.09003  -1.156    0.261    
## Bleomycin:genotypeYFP-ALC1  0.12249    0.08608   1.423    0.170    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1286 on 20 degrees of freedom
## Multiple R-squared:  0.9503, Adjusted R-squared:  0.9429 
## F-statistic: 127.5 on 3 and 20 DF,  p-value: 3.315e-13

cat("AIC: ", AIC(fit2))

## AIC:  -24.70436

simres <- simulateResiduals(fittedModel = fit2)
plot(simres)

fit3 <- lm(NormCounts2 ~ Bleomycin*genotype, data = dataset)
print(summary(fit3))

## 
## Call:
## lm(formula = NormCounts2 ~ Bleomycin * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.15121 -0.03700 -0.01070  0.05473  0.13129 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                 1.00055    0.03746  26.709  < 2e-16 ***
## Bleomycin                  -0.51061    0.03582 -14.255 6.13e-12 ***
## genotypeYFP-ALC1            0.04514    0.05298   0.852    0.404    
## Bleomycin:genotypeYFP-ALC1  0.02064    0.05066   0.408    0.688    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.07569 on 20 degrees of freedom
## Multiple R-squared:  0.9517, Adjusted R-squared:  0.9445 
## F-statistic: 131.5 on 3 and 20 DF,  p-value: 2.479e-13

cat("AIC: ", AIC(fit3))

## AIC:  -50.15771

simres <- simulateResiduals(fittedModel = fit3)
plot(simres)

fit4 <- lmer(Counts ~ Bleomycin*genotype + (1|UID), data = dataset)
print(summary(fit4))

## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ Bleomycin * genotype + (1 | UID)
##    Data: dataset
## 
## REML criterion at convergence: 282.7
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.91400 -0.53151 -0.01178  0.60395  1.83213 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  UID      (Intercept) 55274    235.1   
##  Residual             36717    191.6   
## Number of obs: 24, groups:  UID, 6
## 
## Fixed effects:
##                            Estimate Std. Error       df t value Pr(>|t|)    
## (Intercept)                1271.943    165.583    6.393   7.682 0.000187 ***
## Bleomycin                  -653.663     90.676   16.000  -7.209 2.08e-06 ***
## genotypeYFP-ALC1            417.333    234.169    6.393   1.782 0.121961    
## Bleomycin:genotypeYFP-ALC1 -141.556    128.235   16.000  -1.104 0.285973    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) Blmycn gYFP-A
## Bleomycin   -0.465              
## gntYFP-ALC1 -0.707  0.329       
## Bl:YFP-ALC1  0.329 -0.707 -0.465

cat("AIC: ", AIC(fit4))

## AIC:  294.7363

simres <- simulateResiduals(fittedModel = fit4)
plot(simres)

Quadratic formula

fit5 <- lm(Counts ~ Experiment + poly(Bleomycin,2)*genotype, data = dataset)
print(summary(fit5))

## 
## Call:
## lm(formula = Counts ~ Experiment + poly(Bleomycin, 2) * genotype, 
##     data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -421.90 -106.95  -25.02  149.03  445.60 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            509.21     105.15   4.842  0.00018 ***
## Experimentexp2                         315.50     128.79   2.450  0.02618 *  
## Experimentexp3                         306.88     128.79   2.383  0.02993 *  
## poly(Bleomycin, 2)1                  -1953.50     364.27  -5.363 6.34e-05 ***
## poly(Bleomycin, 2)2                     92.26     364.27   0.253  0.80327    
## genotypeYFP-ALC1                       297.08     105.15   2.825  0.01219 *  
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  -423.04     515.15  -0.821  0.42360    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1  -361.30     515.15  -0.701  0.49316    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 257.6 on 16 degrees of freedom
## Multiple R-squared:  0.8457, Adjusted R-squared:  0.7782 
## F-statistic: 12.53 on 7 and 16 DF,  p-value: 1.973e-05

cat("AIC: ", AIC(fit5))

## AIC:  342.8408

simres <- simulateResiduals(fittedModel = fit5)
plot(simres)

fit6 <- lm(NormCounts ~ poly(Bleomycin,2)*genotype, data = dataset)
print(summary(fit6))

## 
## Call:
## lm(formula = NormCounts ~ poly(Bleomycin, 2) * genotype, data = dataset)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.242931 -0.041198  0.004446  0.065468  0.213157 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           1.000e+00  3.629e-02  27.559 3.58e-16 ***
## poly(Bleomycin, 2)1                  -2.692e+00  1.778e-01 -15.145 1.10e-11 ***
## poly(Bleomycin, 2)2                   1.168e-01  1.778e-01   0.657    0.519    
## genotypeYFP-ALC1                      4.392e-17  5.132e-02   0.000    1.000    
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  3.661e-01  2.514e-01   1.456    0.163    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1 -3.987e-01  2.514e-01  -1.586    0.130    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1257 on 18 degrees of freedom
## Multiple R-squared:  0.9573, Adjusted R-squared:  0.9454 
## F-statistic: 80.71 on 5 and 18 DF,  p-value: 1.097e-11

cat("AIC: ", AIC(fit6))

## AIC:  -24.34165

simres <- simulateResiduals(fittedModel = fit6)
plot(simres)

fit7 <- lm(NormCounts2 ~ poly(Bleomycin,2)*genotype, data = dataset)
print(summary(fit7))

## 
## Call:
## lm(formula = NormCounts2 ~ poly(Bleomycin, 2) * genotype, data = dataset)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.137694 -0.024008  0.002698  0.039744  0.120818 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           0.56680    0.02115  26.794 5.88e-16 ***
## poly(Bleomycin, 2)1                  -1.52597    0.10363 -14.724 1.76e-11 ***
## poly(Bleomycin, 2)2                   0.06622    0.10363   0.639   0.5309    
## genotypeYFP-ALC1                      0.06268    0.02992   2.095   0.0506 .  
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  0.06170    0.14656   0.421   0.6788    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1 -0.24363    0.14656  -1.662   0.1138    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.07328 on 18 degrees of freedom
## Multiple R-squared:  0.9593, Adjusted R-squared:  0.948 
## F-statistic: 84.83 on 5 and 18 DF,  p-value: 7.161e-12

cat("AIC: ", AIC(fit7))

## AIC:  -50.24129

simres <- simulateResiduals(fittedModel = fit7)
plot(simres)

fit8 <- lmer(Counts ~ poly(Bleomycin,2)*genotype + (1|UID), data = dataset)
print(summary(fit8))

## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(Bleomycin, 2) * genotype + (1 | UID)
##    Data: dataset
## 
## REML criterion at convergence: 251.1
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.59469 -0.57627 -0.04532  0.59448  1.51209 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  UID      (Intercept) 54685    233.8   
##  Residual             39074    197.7   
## Number of obs: 24, groups:  UID, 6
## 
## Fixed effects:
##                                      Estimate Std. Error       df t value
## (Intercept)                            716.67     146.57     4.00   4.889
## poly(Bleomycin, 2)1                  -1953.49     279.55    14.00  -6.988
## poly(Bleomycin, 2)2                     92.26     279.55    14.00   0.330
## genotypeYFP-ALC1                       297.08     207.29     4.00   1.433
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  -423.05     395.34    14.00  -1.070
## poly(Bleomycin, 2)2:genotypeYFP-ALC1  -361.30     395.34    14.00  -0.914
##                                      Pr(>|t|)    
## (Intercept)                           0.00811 ** 
## poly(Bleomycin, 2)1                  6.37e-06 ***
## poly(Bleomycin, 2)2                   0.74625    
## genotypeYFP-ALC1                      0.22509    
## poly(Bleomycin, 2)1:genotypeYFP-ALC1  0.30269    
## poly(Bleomycin, 2)2:genotypeYFP-ALC1  0.37624    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) pl(B,2)1 pl(B,2)2 gYFP-A p(B,2)1:
## ply(Blm,2)1  0.000                                  
## ply(Blm,2)2  0.000  0.000                           
## gntYFP-ALC1 -0.707  0.000    0.000                  
## p(B,2)1:YFP  0.000 -0.707    0.000    0.000         
## p(B,2)2:YFP  0.000  0.000   -0.707    0.000  0.000

cat("AIC: ", AIC(fit8))

## AIC:  267.1439

simres <- simulateResiduals(fittedModel = fit8)
plot(simres)

Cubic formula

fit9 <- lm(Counts ~ Experiment + poly(Bleomycin,3)*genotype, data = dataset)
print(summary(fit9))

## 
## Call:
## lm(formula = Counts ~ Experiment + poly(Bleomycin, 3) * genotype, 
##     data = dataset)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -404.54  -95.08   -6.06  104.79  406.96 
## 
## Coefficients:
##                                       Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                            509.208    105.905   4.808 0.000278 ***
## Experimentexp2                         315.500    129.707   2.432 0.029009 *  
## Experimentexp3                         306.875    129.707   2.366 0.032948 *  
## poly(Bleomycin, 3)1                  -1953.495    366.867  -5.325 0.000107 ***
## poly(Bleomycin, 3)2                     92.264    366.867   0.251 0.805088    
## poly(Bleomycin, 3)3                    346.343    366.867   0.944 0.361149    
## genotypeYFP-ALC1                       297.083    105.905   2.805 0.014036 *  
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  -423.045    518.828  -0.815 0.428506    
## poly(Bleomycin, 3)2:genotypeYFP-ALC1  -361.300    518.828  -0.696 0.497592    
## poly(Bleomycin, 3)3:genotypeYFP-ALC1    -1.664    518.828  -0.003 0.997487    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 259.4 on 14 degrees of freedom
## Multiple R-squared:  0.8631, Adjusted R-squared:  0.775 
## F-statistic: 9.804 on 9 and 14 DF,  p-value: 0.0001174

cat("AIC: ", AIC(fit9))

## AIC:  343.9775

simres <- simulateResiduals(fittedModel = fit9)
plot(simres)

fit10 <- lm(NormCounts ~ poly(Bleomycin,3)*genotype, data = dataset)
print(summary(fit10))

## 
## Call:
## lm(formula = NormCounts ~ poly(Bleomycin, 3) * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.10726 -0.04199 -0.01795  0.05358  0.14330 
## 
## Coefficients:
##                                        Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           1.000e+00  2.481e-02  40.307  < 2e-16 ***
## poly(Bleomycin, 3)1                  -2.692e+00  1.215e-01 -22.151 1.97e-13 ***
## poly(Bleomycin, 3)2                   1.168e-01  1.215e-01   0.961 0.350736    
## poly(Bleomycin, 3)3                   4.929e-01  1.215e-01   4.056 0.000918 ***
## genotypeYFP-ALC1                      8.710e-17  3.509e-02   0.000 1.000000    
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  3.661e-01  1.719e-01   2.130 0.049063 *  
## poly(Bleomycin, 3)2:genotypeYFP-ALC1 -3.987e-01  1.719e-01  -2.319 0.033931 *  
## poly(Bleomycin, 3)3:genotypeYFP-ALC1 -1.938e-01  1.719e-01  -1.128 0.276093    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.08594 on 16 degrees of freedom
## Multiple R-squared:  0.9823, Adjusted R-squared:  0.9745 
## F-statistic: 126.5 on 7 and 16 DF,  p-value: 8.384e-13

cat("AIC: ", AIC(fit10))

## AIC:  -41.41727

simres <- simulateResiduals(fittedModel = fit10)
plot(simres)

fit11 <- lm(NormCounts2 ~ poly(Bleomycin,3)*genotype, data = dataset)
print(summary(fit11))

## 
## Call:
## lm(formula = NormCounts2 ~ poly(Bleomycin, 3) * genotype, data = dataset)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.06752 -0.02410 -0.01047  0.03373  0.08122 
## 
## Coefficients:
##                                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                           0.56680    0.01442  39.314  < 2e-16 ***
## poly(Bleomycin, 3)1                  -1.52597    0.07063 -21.605  2.9e-13 ***
## poly(Bleomycin, 3)2                   0.06622    0.07063   0.938  0.36241    
## poly(Bleomycin, 3)3                   0.27939    0.07063   3.956  0.00113 ** 
## genotypeYFP-ALC1                      0.06268    0.02039   3.074  0.00726 ** 
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  0.06170    0.09989   0.618  0.54549    
## poly(Bleomycin, 3)2:genotypeYFP-ALC1 -0.24363    0.09989  -2.439  0.02675 *  
## poly(Bleomycin, 3)3:genotypeYFP-ALC1 -0.09112    0.09989  -0.912  0.37520    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.04994 on 16 degrees of freedom
## Multiple R-squared:  0.9832, Adjusted R-squared:  0.9758 
## F-statistic: 133.7 on 7 and 16 DF,  p-value: 5.446e-13

cat("AIC: ", AIC(fit11))

## AIC:  -67.4721

simres <- simulateResiduals(fittedModel = fit11)
plot(simres)

fit12 <- lmer(Counts ~ poly(Bleomycin,3)*genotype + (1|UID), data = dataset)
print(summary(fit12))

## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Counts ~ poly(Bleomycin, 3) * genotype + (1 | UID)
##    Data: dataset
## 
## REML criterion at convergence: 221.8
## 
## Scaled residuals: 
##      Min       1Q   Median       3Q      Max 
## -1.55170 -0.46246 -0.07901  0.43109  1.35742 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  UID      (Intercept) 55544    235.7   
##  Residual             35638    188.8   
## Number of obs: 24, groups:  UID, 6
## 
## Fixed effects:
##                                       Estimate Std. Error        df t value
## (Intercept)                            716.667    146.575     4.000   4.889
## poly(Bleomycin, 3)1                  -1953.495    266.975    12.000  -7.317
## poly(Bleomycin, 3)2                     92.264    266.975    12.000   0.346
## poly(Bleomycin, 3)3                    346.343    266.975    12.000   1.297
## genotypeYFP-ALC1                       297.083    207.289     4.000   1.433
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  -423.045    377.560    12.000  -1.120
## poly(Bleomycin, 3)2:genotypeYFP-ALC1  -361.300    377.560    12.000  -0.957
## poly(Bleomycin, 3)3:genotypeYFP-ALC1    -1.664    377.560    12.000  -0.004
##                                      Pr(>|t|)    
## (Intercept)                           0.00811 ** 
## poly(Bleomycin, 3)1                  9.26e-06 ***
## poly(Bleomycin, 3)2                   0.73563    
## poly(Bleomycin, 3)3                   0.21892    
## genotypeYFP-ALC1                      0.22509    
## poly(Bleomycin, 3)1:genotypeYFP-ALC1  0.28444    
## poly(Bleomycin, 3)2:genotypeYFP-ALC1  0.35748    
## poly(Bleomycin, 3)3:genotypeYFP-ALC1  0.99656    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##             (Intr) pl(B,3)1 pl(B,3)2 pl(B,3)3 gYFP-A p(B,3)1: p(B,3)2:
## ply(Blm,3)1  0.000                                                    
## ply(Blm,3)2  0.000  0.000                                             
## ply(Blm,3)3  0.000  0.000    0.000                                    
## gntYFP-ALC1 -0.707  0.000    0.000    0.000                           
## p(B,3)1:YFP  0.000 -0.707    0.000    0.000    0.000                  
## p(B,3)2:YFP  0.000  0.000   -0.707    0.000    0.000  0.000           
## p(B,3)3:YFP  0.000  0.000    0.000   -0.707    0.000  0.000    0.000

cat("AIC: ", AIC(fit12))

## AIC:  241.8309

simres <- simulateResiduals(fittedModel = fit12)
plot(simres)

Compare Results

ICtab(fit1,fit2,fit3,fit4,
      fit5,fit6,fit7,fit8,
      fit9,fit10,fit11,fit12,
      base=T)

##       AIC   dAIC  df
## fit11 -67.5   0.0 9 
## fit7  -50.2  17.2 7 
## fit3  -50.2  17.3 5 
## fit10 -41.4  26.1 9 
## fit2  -24.7  42.8 5 
## fit6  -24.3  43.1 7 
## fit12 241.8 309.3 10
## fit8  267.1 334.6 8 
## fit4  294.7 362.2 6 
## fit1  339.7 407.2 7 
## fit5  342.8 410.3 9 
## fit9  344.0 411.4 11

Final Result

fit <- fit11

output <- coef(summary(fit))
output <- output[grep("Bleomycin", rownames(output)),]

rownames(output) <- gsub("poly\\(|, [1-3]\\)","", rownames(output) )
rownames(output) <- gsub("genotype",  paste0(" ",levels(dataset$genotype)[1], " vs. "), rownames(output))
rownames(output)[!(grepl("vs", rownames(output)))] <- paste(rownames(output)[!(grepl("vs", rownames(output)))], levels(dataset$genotype)[1],  sep = " in " )


# suggested result table
kable(output, row.names = T)

	Estimate	Std. Error	t value	Pr(>|t|)
Bleomycin1 in WT	-1.5259664	0.0706299	-21.6050932	0.0000000
Bleomycin2 in WT	0.0662200	0.0706299	0.9375630	0.3624070
Bleomycin3 in WT	0.2793890	0.0706299	3.9556733	0.0011329
Bleomycin1: WT vs. YFP-ALC1	0.0616959	0.0998858	0.6176642	0.5454882
Bleomycin2: WT vs. YFP-ALC1	-0.2436254	0.0998858	-2.4390393	0.0267529
Bleomycin3: WT vs. YFP-ALC1	-0.0911161	0.0998858	-0.9122025	0.3752023

write.table(output, file = "FigureS5L_Stats.txt", quote = F, sep = "\t", row.names = T, col.names = NA)

Anova

fit11a <- lm(NormCounts2 ~ poly(Bleomycin,3)*genotype, data = dataset)
fit11b <- lm(NormCounts2 ~ poly(Bleomycin,3)+genotype, data = dataset)

# anova table
anova(fit11a, fit11b)

## Analysis of Variance Table
## 
## Model 1: NormCounts2 ~ poly(Bleomycin, 3) * genotype
## Model 2: NormCounts2 ~ poly(Bleomycin, 3) + genotype
##   Res.Df      RSS Df Sum of Sq      F Pr(>F)
## 1     16 0.039909                           
## 2     19 0.057774 -3 -0.017865 2.3875 0.1071